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Abstract 

We analyze canonical partition sums with a fixed quark number q in 
lattice QCD. In particular the fermion determinant is projected to consist 
of loops with a net winding of q. We compare different strategies for the 
evaluation of these canonical determinants and discuss their properties 
using pure SU(3) gauge configurations for a wide range of temperatures. It 
is demonstrated that the center sectors of the Polyakov loop very strongly 
manifest themselves in the behavior of the canonical determinants in the 
deconfined phase, and we discuss physical implications of this finding. 
Furthermore we show that canonical determinants with even numbers of 
quarks become enhanced in the high temperature phase. 



1. Introductory remarks 



With the already running and planned experiments at RHIC, LHC and FAIR, 
Quantum Chromodynamics with temperature and non-vanishing density has 
recently shifted into the center of attention. Of particular interest are non- 
perturbative results focusing on various aspects, such as obtaining information 
about the phase diagram, characterization of the different phases, and a possible 
understanding of the mechanisms that drive the various transitions. In principle 
lattice QCD is such a non-perturbative approach, but numerical simulations at 
finite density are plagued by serious phase cancellation problems. 

There exist two alternative approaches to QCD with finite density: One can 
work grand-canonically using a chemical potential, or alternatively consider the 
canonical ensemble, i.e., a fixed number of quarks. The majority of lattice calcu- 
lations at finite density uses a chemical potential, while the canonical approach 
has seen comparatively little attention [IJ - [6J. 

in this article we revisit the canonical approach. Using a recently proposed 
dimensional reduction for the fermion determinant [6] we can numerically eval- 
uate canonical determinants much faster, i.e., fermion determinants describing 
a fixed number of quarks. We combine the dimensional reduction with different 
strategies for computing canonical determinants and compare their efficiency. 

Canonical determinants are not only objects that might provide new and 
more efficient strategies for lattice QCD with finite density, but they also have 
interesting properties themselves. In particular they show simple transformation 
properties under center transformations - a symmetry that is of crucial impor- 
tance for understanding the high temperature phase transition of pure gauge 
theory. Thus the numerical analysis of canonical determinants might provide 
important insight into mechanisms responsible for transitions in the QCD phase 
diagram. In this article we demonstrate that the influence of center symmetry 
and its breaking are very strongly manifest in the properties of the canonical 
determinants in the high temperature phase and we discuss the physical im- 
plications of this finding. Furthermore, we show that in the high temperature 
phase the canonical determinants with odd quark numbers become suppressed. 

2. Canonical partition sums on the lattice 

In this section we present our lattice Dirac operator and compile formulas for 
canonical partition sums and canonical determinants. This serves to set our 
notational conventions and to summarize the formulas which we later need for 
the analysis of the canonical determinants. 
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2.1. Wilson fermions with chemical potential 



In this article we use Wilson's lattice Dirac operator. The chemical potential is 
introduced as a boundary condition, such that the lattice Dirac operator reads 

DinUy = 5,,y - Ui{x) 5,^3.^ (1) 

i=±i ^ 

n=-P/2+l 

13/2 ^^^^ 

— 2 ^4(iC) 1)^ ^a;4,n ^j/4,n-l 

n=-/3/2+2 

+ e^"'^ K U^ix, Sg^ff '5x4-/3/2+1 ^y4,P/2 ■ 

The lattice indices x = {x,Xi),y = {y,y4) run over the x (3 lattice, by j 
we denote the shift vector in j'-direction, and we use the convention U-j{x) = 
Uj{x — j)^ ■ P is the inverse temperature in lattice units, i.e., /? is a positive 
integer. The temporal lattice index X4 runs from — /3/2 + 1 to /?/2. In our 
notation the hopping parameter k is related to the bare quark mass m via k = 
1/(4 + m). The chemical potential n is attached to the link terms connecting 
■^4 = /3/2 and .T4 = —(3/2 + 1 in the usual exponential form. These terms 
also come with an additional relative sign which implements the anti-periodic 
temporal boundary conditions for the fermions. The link variables Uy{x) are 
elements of the gauge group SU(3). The corresponding color indices, as well as 
the Dirac indices are not displayed explicitly and we use vector/matrix notation 
for them. 

For later use we recollect two important properties of the lattice Dirac op- 
erator -D(/u) and its determinant: The first one is its behavior under a similarity 
transformation with 75 (^u is real here), 

75i^(/x)75 = (2) 

which often is referred to as generalized 75-hermiticity. It implies for the deter- 
minant of the Dirac operator the identity 

det[Z)(^)] = det[L»(-//)]* , (3) 

where the asterisk denotes complex conjugation. 
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The second property reflects the behavior under a time reflection of the 
fields, and is conveniently expressed by the identity 



det 



ru 



det[Z)(-/i)] 



(4) 



On the left-hand side the determinant is evaluated for the time reflected gauge 
field TU (see, e.g., [7j for a detailed definition of the time reflection), while on 
the right-hand side the original gauge field configuration U is used. Combining 
([3]) and (jlj) gives rise to 



det[D(^)] 



TU 



det[Z)(^)]* 



(5) 



2.2. Definition of canonical partition sums 

This subsection gives a short summary of the definition of canonical partition 
sums to set our notation. For N mass degenerate flavors of fermions, the grand 
canonical partition sum Zccif^) 3t chemical potential /i is given by 

ZgM = (^{det[D{f,)]f^^ . (6) 

In the following we will often refer to the fermion determinant det[L'(/u)] or its 
A^-th power as the grand canonical determinant. By {X)c we denote the path 
integral expression for some observable X in pure gauge theory, 

{X)g = J V[U]e-'olu]x[u] , (7) 

where Sg[U] is the gauge action and ^[11] denotes the path integral measure 
over all gauge fields, which on the lattice is given by the product of the Haar 
measures for all link variables Ui,{x). 

For the grand canonical partition sum one may perform a fugacity expan- 
sion, i.e., write it as a sum over canonical partition sums ^(9) with a fixed 
quark number q £ Ij, 

ZgM = Y.^'^'^Z^"'^ ■ (8) 

As a matter of fact, for a finite x f3 lattice the sum in the fugacity expansion 
runs only over the finite interval q G [—qmax,qmax], where, e.g., for one flavor 
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of Wilson fermions one has qmax = 6L^. The canonical partition sums Z^''^ 
with a fixed quark number q may themselves be written as path integrals, 

= (det^[Z)]W)^, (9) 

where the canonical fermion determinant with a fixed quark number q is de- 
noted as det AT [D] . We attach the number of flavors as a lower index 
and omit the argument fi of the Dirac operator, since it has disappeared from 
the canonical expressions and enters only through the expansion coefficients 
exp(/xg/3) in the fugacity expansion ([8]). 

We stress at this point, that the canonical A^-flavor determinant det^lD]^'^^ 
with fixed quark number q does not factorize into one-flavor canonical determi- 
nants. This is different from the grand canonical determinant which, as seen in, 
e g-, dS]), is simply obtained by taking the A^-th power of the grand canonical 
one-flavor determinant det[D{fi)]. 

We remark that once an efficient way of computing the fermion determi- 
nant detjv ("^^ for a fixed quark number q is available, the expression for the 
corresponding canonical partition sum Z^^^ is rather simple. Eq. ([9]) shows that 
up to a normalization it can be written as an expectation value (...)g in Pure 
gauge theory. 

A possible way of computing detN[D]^'^^ is through a Fourier transformation 
with respect to imaginary chemical potential, 

detNlD]^'^ = ^ £dipe'''^'^ (det[Z)(/x = . (10) 

The explicit form ([T]) of the lattice Dirac operator shows that also imaginary 
chemical potential is just a boundary condition which multiplies the forward 
running temporal link variables of the last time slice with the phase e*"^, while 
the backward running temporal link variables on the last time slice are multiplied 
with e-^'f. 

Combining Eqs. ([6j) - (ITOll . one obtains the fugacity expansion 

[det[D{f,)]y = ^ e^^"^ detNlD]^'^ , (11) 
1 

of the grand canonical determinant {det[D{fi)])^ in terms of the canonical 
determinants detN[D]^'^\ 

Equation (fTTIl may now be used to obtain an expression of the canonical 
determinants with > 1 flavors in terms of the canonical one-flavor determi- 
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nants. Exploring ()TT|1 for = 1 one finds 

/ \ N 



(det[D(/x)])'^ = (^e'^^^detiiD]^^)) (12) 



N 

<?i 92 gjv i=i 

N 

= E e^'"" E - E '^fe + - + Q) n deti[I)](''^), 
q <?i giv i=i 

where 6{qi + ... + qN,Q) denotes a Kronecker delta which enforces qi + q2 + 
... + qN = q. Comparison with (fTTi) allows one to read off the expression 

AT 

det^[D](^) = ^5]...5]%i + ... + g^,g)ndeti[Z)](*), (13) 

91 12 QN j=l 

which we will later use to compute the canonical determinants detjv[-D](^) for 
the case of > 1 flavors from the one-flavor canonical determinants deti[L']'''^^. 
We conclude this subsection with exploring the implications of 75-hermiticity 
and the time reflection symmetry ([4]) for the canonical determinants: Com- 
bining 75-hermiticity in the form of ([3]) with the fugacity expansion (jlll) . one 
finds the relation 

det N[D]^''i^ = (det N[D]^'^^y , (14) 

between canonical determinants with positive and negative quark numbers q. 
For vanishing quark number q = the canonical determinant is real, i.e., 

det^[D](0) G R. 

In a similar way we use ([4]) together with the fugacity expansion (jlll) and 
arrive at 

det^[D](^)|^^ = {detNlD]^^^ \J , (15) 

which shows that under a time reflection of the gauge configuration canonical 
determinants turn into their complex conjugate. Since in the path integral {■..)g 
over all gauge fields both, a configuration U and its time reflected image TU 
appear with the same weight, we find 

Z^i^ = (det^[D](^)) = ^(detjv[-D]^^) + det jv [I?] 



u 



rui G 



RedetNiDf^') G R, (16) 
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i.e., the canonical partition sums are real. Via the fugacity expansion ([8]) 
this implies that also the grand canonical partition sum Zqc is real. 

2.3. The fermion loop picture for canonical determinants 

The canonical fermion determinants detAr[D]('^) have an interesting interpreta- 
tion as collections of fermion loops with a fixed total winding number q: As a 
gauge invariant quantity on the lattice, the grand canonical fermion determinant 
(det[Z)])^ can be decomposed into sums over products of closed loops on the 
lattice dressed with the gauge variables on the links the loops consist of. One 
finds 

{dei[D]f = ^c,n (Tr I] U,{x)\ . (17) 

In our notation s runs over sets S of closed loops, Cg is a complex valued 
coefficient, and / runs over the individual loops in the set s. The second product 
is over the links (x,/i) the loop I is made of, and the trace is over color indices. 
One may, e.g., use hopping expansion to compute some terms in (fT7|) explicitly. 

Since the loops I are closed, one can determine their winding number qi 
around compactified time with qi £ (negative values correspond to backward 
winding). An individual term in the sum (fT7i) . which consists of a set s of loops, 
may then be assigned a total winding number q^ = J2i£s H- 

When one uses the definition (jlOp of the canonical determinants detN[D]^'^\ 
a boundary condition phase e^"^^ is attached to the temporal link variables on 
the last time slice. Thus a loop / that winds qi times picks up a factor of e*"?'*^. 
An individual term in the sum (fTTi) which is made from a product of loops in 
a set s thus acquires a factor of e*'^"'^. The Fourier integral in (flO]l projects to 
a fixed total winding number qg = q, such that for the canonical determinants 
we obtain the loop representation 

detr,[D]^^^ = E [Tr n U,{x)] , (18) 

where S^'^^ now contains only sets s of loops that have a total winding number 
qs = q. Thus the canonical determinant det at [-D] ^'^^ is made of products of 
loops with a total winding number of q. 
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2.4. Center transformation properties of the canonical determinants 

SU(N) gauge theory has an interesting property which plays an important role 
for the QCD phase transitions, namely a symmetry under center transformations. 
In lattice language the center transformation consists of multiplying all temporal 
link variables at a fixed time slice X4 = x| with a fixed element z of the center 
of the gauge group, 

C/4(x, X4) > X4), (19) 

while all other link variables remain unchanged. For the case considered here, 
i.e., gauge group SU(3), the center is Z3 and thus z has the three possible 
values 

z G {1 , e'2-/3 , e-*2V3| _ (20) 

The gauge action Sg^\ is invariant under such a transformation, not, however, 
the grand canonical fermion determinant. 

According to Eq. p8]l the canonical determinants det at [D] '^'^^ can be written 
as sums over products of loops, where each term has a total combined winding 
number of q for the loops the term consists of. Since the transformation (fl9]) 
affects only a single time slice, the total factor that is picked up by the loops is 
given by z"?. Thus we conclude that under the center transformation the 
canonical determinants transform as 

det7v[I)]^^^ ^ z^detTvl^*]^"^ = z^^o^^^etjvi-Dl^^) , (21) 

i.e., they pick up a phase which depends on the quark number q. Since z may 
only assume the values listed in ([^U|) , which form the center group Z3, the 
powers z'^ are also restricted to these values. In the last expression of (l2TIl we 
have stressed this fact by introducing the concept of triality, defined as q mod 3. 
Equation (|2T1) shows that only the zero triality determinants remain invariant, 
while for non-zero triality a phase is picked up under the center transformation. 

The center symmetry has a simple physical consequence: Performing a cen- 
ter transformation of the integration variables U in the path integral ([7]) leaves 
both the gauge action and the integration measure invariant. Thus as long as 
the center symmetry is not broken spontaneously in {...)g, the only non-trivial 
transformation of an expectation value can come from the integrand inserted in 
{...)g- For the canonical partition sums this leads to 

= (det^l^l^"))^ = z«(det^[D](«))^, (22) 

which implies that Z^'^^ must vanish unless z'^ = 1. Thus only if the canonical 
partition sum Z^'^^ is of zero triality, i.e., gmodS = 0, it is non-vanishing. 
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To include this finding in the fugacity expansion ([8]), we write the canonical 
partition sum as 

Zgc = e^'^^I^Z^^^^ , (23) 

b 

where we take into account only the non-vanishing zero triality terms with q = 3b 
by introducing the baryon number 6 € We stress again, that this chain of 
arguments holds only if the center symmetry is not broken spontaneously by 
{...)g- The impact of a spontaneous breaking will be discussed in Section 5. 

3. Techniques for computing canonical fermion determinants 

In the last section we have addressed two possible ways of computing the canon- 
ical determinants det]\f[D]^'^\ Equation (jlOp defines them through a Fourier 
transformation, while (jl8l) formally expresses them as a sum over loops with a 
fixed total winding number. Both equations are not very suitable for an immedi- 
ate application in a numerical evaluation, and in this section we discuss various 
strategies for an efficient evaluation of canonical determinants. We begin this 
discussion with a presentation of a dimensional reduction formula [6] for the de- 
terminant that will be a powerful tool for all approaches which we subsequently 
discuss. We will restrict ourselves to the simple case of only one flavor, since 
the case of more than one flavor, i.e., > 1, can be reduced to = 1 with 
the help of (fTSll . 

3.1. Dimensional reduction of the fermion determinant 

For later use we summarize the essential formulas for an exact dimensional 
reduction of the fermion determinant [6]. For that purpose the lattice is de- 
composed into 4 domains, A^'-\i = 1,2,3,4, which are distinguished by the 
values of the temporal lattice index X4. In particular A*^^) consists of all lattice 
points in the range from X4 = — /3/2 -|- 1 to — 1 and A^^) has X4 running from 
-|-1 to P/2 — 1. The domains A^^^ and A^^^ consist of single time slices with 
X4 = and X4 = P/2, respectively. The domain decomposition of the lattice is 
illustrated in Fig. [H 

Based on the domain decomposition the Dirac operator is split into several 
pieces: The terms D^^\i = 1,2,3,4, live on the individual domains A^^\i = 
1,2,3,4. The terms D^^'^), d(2.3)^ £)(3,4) j-)i4,i) contain the forward hop- 
ping temporal terms that connect the domains, and D^'^'^\ D^^''^\ £)(4,3) 
D^^''^^ are the corresponding backward hopping terms. 
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Figure 1: Partition of the lattice into four domains A^'^\i = 1,2,3,4: The 
vertical direction in the plot is time, while a line of dots in horizontal di- 
rection represents a whole time slice. On the left-hand side we indicate the 
values of the time argument. Lattice sites in a common domain are enclosed 
in a dashed contour, and we label the domains on the right-hand side. 

The fermion determinant may be factorized exactly [6] in the form 

det[D{n)] = Ao det [1 - Hq - e^^H+i - e'^'^H.i] , (24) 

where the terms Hq,H±i live on only the single time slice that constitutes the 
domain A^"^). They are made out of the Dirac operator terms D^^\D^'^'^^ and 
the corresponding propagators, 

= 5(4)5(^-2)5(2)5(2'^), (25) 

where 

5(1) = ^ ^(3) ^ p(3))-l 

^(2) ^ ^(2) _ ^(2,1)^(1)^(1,2) _ ^(2,3)^(3)^(3,2)^ 

^(4) = Z){4)_Z)(4.1)S(1)Z)(1.4)_ ^(4,3)^(3)^(3,4)^ 

5(2) = (5(2))-i , 5(4) = (5(4))-i , 
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= d(^'^)5(i)z)(1'2) , 

^(4,2) ^ ^(4,3)^(3)^(3,2) ^ 
^(2,4) ^ ^(2,1)^(1)^(1.4) 

^(2,4) ^ ^(2,3)^(3)^(3,4) ^ (26) 

The subscript attached to the terms D^^'^^ and D^^'^^ indicates that there 
the chemical potential was set to 0, and instead is taken onto account by the 
exponential factors attached to the terms H±i in Eq. (j24l) . The overall factor 
Aq in (I24p is independent of the chemical potential fi and given by 

^0 = det [Z)(i)] det [d'^^^] det [5(2)] det [5^^)] . (27) 

We stress again that the formula ([24]) is an exact factorization that does not 
contain any approximation and holds for arbitrary background gauge fields. It 
provides a powerful tool for the analysis of lattice QCD with chemical potential: 

• Eq. (l24ll factorizes the piece that depends on the chemical potential fi 
and reduces it to a determinant which lives on a single time slice. Due to 
this dimensional reduction the numerical evaluation of the ^-dependent 
piece is considerably less demanding in a numerical analysis. 

• The ^-independent and real factor Aq was shown [6J to carry most of the 
scale, i.e., the overall size of the determinant which rapidly grows with 
the volume and shows a strong quark mass dependence. The second, 
//-dependent determinant in (l24ll remains small compared to Aq, and in 
isolated form contains the terms that are problematic in lattice QCD with 
chemical potential, i.e., the complex phase. 

• Finally the factorization p4|) allows one to compute the canonical deter- 
minants using a generalized hopping expansion with excellent convergence 
properties. We will elaborate on this approach below. 

3.2. Fourier projection method 

The most straightforward method for the determination of the canonical deter- 
minants is the Fourier projection as shown in Eq. (flOl) . As a matter of fact 
it is sufficient to apply this formula to the one-flavor case = 1, since the 
determinants for > 1 may be built up using ([13]). 

A problem of the Fourier projection method is that for larger values of q the 
Fourier integral has to be computed very precisely. For the numerical evaluation 
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this implies that the determinant has to be computed for many values of 99 in the 
interval [—it, it]. Since the numerical cost for the evaluation of a determinant 
increases with the size of an iV x matrix like A^^, this is a rather expensive 
approach. 

The situation improves considerably when the factorization formula dis- 
cussed in the last subsection is applied [6J. Using (j24l) together with (jlOl) 
one finds 



In this expression the determinant is for matrices that live on only a single time 
slice. Thus the size N of the matrix is reduced by a factor of and the cost 
for computing the determinant goes down by a factor of 1//?^. For lattice sizes 
typically used in finite-/^ calculations one has /3 = 4 - 6, such that the speed 
up of the calculation is of order 0(100). Thus, at a fixed numerical cost a 
much larger number of intermediate values of ip becomes available and Fourier 
projection becomes an attractive approach also for larger values of q. 

3.3. Generalized hopping expansion 

It may be shown that the terms Hq,H±i consist of paths of at least j3 links, 
dressed with link variables [6J. Consequently a product of n terms Hq,H±i is 
made of paths with a length of at least n(3 links. Since products of link variables 
along such paths decrease in amplitude roughly exponentially with the length 
of the path, this property can be used to set up an efficient expansion [6J of 
the determinant, which may be viewed as a generalization of the usual hopping 
expansion. Using det[M] = exp(TrlnM) one finds 



deti 



rdif e-^'^'P det [1 - - e'^PR+i - e-'^H_i] . (28) 



det[D{iJ,)] 



Aq det [1-Ho- ef'^H+i - e'^'^H^i] 



(29) 




where we have introduced 



00 -| 



(30) 



n=l ki+...+k„=q 
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For a perturbative expansion the sum in the exponent of Eq. ([29|) has to be 
cut off such that q is restricted to g e [—Q:Q] (see last expression in (I29p ). 
The exponential in the last line of ([29ll may then be expanded and the resulting 
series can be organized with respect to powers of exp(ib^/3). Comparison with 
(jlip for the case of iV = 1 allows one to read off the canonical determinants 
deti[D]('^) as the coefficients of the terms exp(^q/3). In a practical application 
also the infinite sum (l30]l has to be truncated at some finite value of n. In [6j it 
was demonstrated that the contributions decrease exponentially with increasing 
n, such that the series (l30]l converges quickly and the effect of the truncation 
can be controlled systematically. One finds that the results for deti[D]'^'^^ from 
the expansion are reliable up to \q\ < Q (see also Section 4). 

3.4. Generalized Fourier method 

Yet another method for the evaluation of the canonical determinants was sug- 
gested in [2]. The basic idea is to expand the logarithm of the determinant into 
its Fourier components and then to expand the exponential of this series. In our 
presentation here we directly apply the method to the dimensionally reduced 
determinant (I24p . Evaluating this expression at imaginary chemical potential, 



det[D{fi = iip/13)] = Aa det [1 - - e'^^H+i - e-''^H_i] , (31) 



we define the expansion coefficients c„, where cq is real, while the c„ with n > 
are complex and the asterisk again denotes complex conjugation. These coef- 
ficients may, e.g., be computed by a numerical evaluation of the corresponding 
Fourier transforms of lndet[l - Hq - e*^F+i - e"*'^i?_i]. 

An alternative strategy for the determination of the c„ is provided by the 
generalized hopping expansion discussed in the last section. A comparison of 
(|3T]1 with the last expression of ([29]1 shows that the coefficients c„ can be 
identified as 





(33) 





(32) 




Q 
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where we have cut off the sum over all Fourier coefficients at some maximal value 
Q. Describing the coefficients by their absolute value and phase, introduced as 

Cn = with dn = 2\Cn\ , (34) 

we find upon inserting this form in (|33p 

det^f^l^") ~ ^ j^d^e-^'^'P exp i^o + cos(nv; + 5„)j 

oo oo °° ( ^ Id \^^\ 

fc2=0 fc3=0 kM=0 \n=2 / 

•'-^ \n=2 

Writing the cosines in the last factors again as exponentials, one can solve the 
remaining (/9-integrals explicitly using 

^ ^gdicos{^+5i)g-im<^ ^ ^imS^l^^di), (36) 
-vr 27r 

where the Im denote the Bessel functions of the first kind. As for the generalized 
hopping expansion, the necessary truncations used in the generalized Fourier 
method may be controlled systematically due to the decrease of the coefficients 
Cn (i.e., the — T^")) with increasing n. 



4. Numerical evaluation of canonical fermion determinants 

Having presented various strategies for the numerical evaluation of the canonical 
determinants in the last section, we now briefly discuss our experience with 
the different approaches. The tests were done for the free case, as well as 
for SU(3) configurations from the pure gauge ensemble ([7]). We begin with 
briefly describing these gauge ensembles in the next subsection and subsequently 
present the outcome of our comparison. 



4.1. The gauge ensembles 

It was already stressed in Section 2.2 that in the canonical approach the canon- 
ical partition sums Z'^''^ may be expressed ([9|) using the expectation value (...)g 
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of pure gauge theory as given in ([7]). In this paper we perform various numerical 
tests using such pure gauge theory ensembles. These tests are done on three 
different lattice volumes, 6^ x 4, 8^ x 4 and 10^ x 4. For the former two we 
use ensemble sizes of 500 (1000 for two ensembles) configurations, while for the 
largest lattice 100 configurations were evaluated. We work with the Liischer- 
Weisz gauge action [S] and for an estimate of the scale we in some places use 
the lattice spacing determined in [9] from the Sommer parameter. The update 
was implemented using a mix of overrelaxation and Metropolis sweeps. Our up- 
date furthermore contains a random center rotation of the configurations. As 
discussed, this is always a symmetry of the gauge theory as long as the system 
size is finite. Since we are interested in the behavior of the canonical determi- 
nants in all three Polyakov loop sectors (compare Section 2.4), we decided to 
include a random center rotation in the update in order to generate an equal 
share of configurations in each sector. 

For the two smaller lattice sizes we use several different values of the gauge 
coupling /3g to scan the low and high temperature regimes. For our lattices 
with time extent /? = 4, the critical value of the gauge coupling is near (3g = 7.8. 
We use the analysis of the critical temperature for the LUscher Weisz gauge 
action in pTO] to determine the ratio T/Tc for our ensembles, which we list in 
Table 1, together with other properties of the gauge configurations. 

For later use and in order to illustrate properties of our configurations, in 
Fig. [2]we show a distribution of the values of the Polyakov loop in the complex 
plane. It is obvious that for low temperature (i.e., the smaller values of the gauge 
coupling (5g) the Polyakov loop values scatter isotropically near the origin in the 
complex plane. For high temperature the well known anisotropical distribution 
starts to form: The data points move away from the origin and the phases Op 
of the Polyakov loop are close to the three center values, i.e.. Op ~ O,0p ~ 
±27r/3. As long as the lattice is finite, the Polyakov loop sectors corresponding 
to these three values are completely equivalent. In the infinite volume limit, 
the underlying center symmetry (see Section 2.4) is broken spontaneously in 
pure gauge theory at high temperature and the system selects one of the three 
sectors. 

In most cases the quark mass m in the Dirac operator was set such that it 
corresponds to 100 MeV on the scale of the pure gauge theory. For two of our 
ensembles we also studied additional quark masses of 50 and 200 MeV. The 
parameters of our gauge ensembles are summarized in Table [TJ 

The determinant calculations and matrix inversions necessary for computing 

^We use the subscript G to distinguish the inverse gauge coupling /3g from the inverse 
temperature f3. 
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Figure 2: Scatter plots of Polyakov loop values in the complex plane for 500 
configurations on our 8^ x 4 lattices at various values of the inverse gauge 
coupling Pg- 



the fermion determinant according to Equations (l24l) - (l27ll are implemented 
with standard LAPACK and BLAS routines. All errors we show are statistical 
errors determined with single elimination Jackknife. 

Concluding this section we need to stress an important caveat that must 
be kept in mind when using the pure gauge theory configurations. Although 
formally one obtains the correct canonical partition function from the expecta- 
tion value, Z^'^^ = {det i\f[D]^'^^)G, a numerical approach might be cumbersome. 
In particular one may expect that a Monte Carlo simulation which uses only 
the weight exp(— S'g'[[/]) of the pure gauge theory will predominantly produce 
configurations with rather small values of det jv [^] '-'^^ and thus only a few 
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size 




rr 1 
a [rmj 


1 [McV J 




m [McV J 


if conis. 


6x4 


7.4 


0.234 


210 


0.70 


100 


500 


6x4 


7.5 


0.213 


232 


0.77 


100 


500 


6x4 


7.6 


0.194 


254 


0.85 


100 


500 


6x4 


7.7 


0.177 


279 


0.93 


100 


500 


6x4 


7.8 


0.161 


306 


1.02 


100 


500 


6x4 


7.9 


0.148 


334 


1.11 


100 


500 


6x4 


8.0 


0.135 


364 


1.21 


100 


500 


6x4 


8.1 


O.lzo 


396 


1.32 


100 


500 


6x4 


8.2 


0.115 


430 


1.43 


100 


500 


8x4 


7.4 


0.234 


210 


0.70 


50, 100, 200 


1000 


8x4 


7.5 


0.213 


232 


0.77 


100 


500 


8x4 


7.6 


0.194 


254 


0.85 


100 


500 


8x4 


7.7 


0.177 


279 


0.93 


100 


500 


8x4 


7.8 


0.161 


306 


1.02 


100 


500 


8x4 


7.9 


0.148 


334 


1.11 


100 


500 


8x4 


8.0 


0.135 


364 


1.21 


100 


500 


8x4 


8.1 


0.125 


396 


1.32 


100 


500 


8x4 


8.2 


0.115 


430 


1.43 


50, 100, 200 


1000 


10 X 4 


7.5 


0.213 


232 


0.77 


100 


100 


10 X 4 


8.1 


0.125 


396 


1.32 


100 


100 



Table 1: Table with the parameters of our gauge configurations. 

configurations could dominate the result. The situation might be improved 
by implementing, e.g., a reweighting approach in a fully dynamical simulation. 
However, already the approximation allows for interesting insights into the be- 
havior of the canonical determinants and we expect that in particular the center 
symmetry properties will carry over to simulations with more advanced sampling 
techniques. 

4.2. Comparison of the numerical approaches to canonical determinants 

In order to assess the three different strategies for computing the canonical 
determinants which we discussed in Sections 3.2 - 3.4, we compared the three 
methods in the free case and for pure gauge theory ensembles on the 6^ x 4 and 
8^ x4 lattices. In all three approaches we made use of the dimensional reduction 
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formula ([24l) . Since this is an exact expression which may be evaluated very 
precisely |6J, this step does not introduce any approximation, but considerably 
speeds up all three approaches. 

Using the standard Fourier method of Section 3.2 one can in principle eval- 
uate the canonical determinants for arbitrary large q, as long as the numerical 
integration scheme used for the Fourier transformation is accurate enough. In 
our implementation we experimented with the trapezoidal integration scheme, 
as well as Simpson's rule, and found almost no difference for the cases we stud- 
ied. Using a sufficiently large number N^p of intermediate values of € [— vr, vr], 
one may compute the canonical determinants up to g = 10 with almost arbi- 
trary precision (on 8^ x 4). We considered N^p up to N^p = 4096, which for 
the smaller values of q is certainly a very conservative number. On the other 
hand, for the larger values of q such an order of magnitude for A'^^ is required 
for sufficient accuracy. By increasing the number A'^^ of points used for the in- 
tegration one can systematically improve the accuracy for each volume and set 
of parameters ficrn, q. Thus we consider the results from the Fourier approach 
as our standard values where we compare the other results to. 

An independent check of the accuracy for the numerical evaluation of the 
canonical determinants is obtained by summing the fugacity expansion (jlll) . As 
we will discuss in Section 5, for the free case and also for our ensembles at 
the lower values of /3g, we found an excellent agreement between the sum and 
the exact expression for the grand canonical determinant. For larger (ic higher 
orders of q are needed (compare Fig. [8]) which were not calculated in this work. 

The other two methods, the generalized hopping expansion (Section 3.3) 
and the generalized Fourier method (Section 3.4) have a cutoff parameter which 
plays the same role in both approaches. This is the number Q in Equations (l29ll 
and (i33]l . respectively. The parameter Q sets the cutoff up to which winding 
number the terms in the logarithm of the determinant are taken into account. 
The actual calculation of the terms is, however, different in the two approaches: 
In the generalized hopping expansion they are obtained as sums over traces of 
products of Hq,H±i, while for the generalized Fourier method they are again 
obtained through a numerical Fourier integral. 

In order to obtain the canonical determinants from the expansion of the 
logarithm of the grand canonical determinant, also the exponential has to be 
expanded. For the two methods described in Sections 3.3 and 3.4 this second 
expansion is done in different ways, but in principle, both methods allow to 
formally calculate canonical determinants det at [D] for arbitrary large q in the 
second expansion, irrespective of the cutoff value Q. However, in our compar- 
ison we found that for both methods sizable deviations from the results of the 
standard Fourier method appear as one increases q beyond the cutoff Q. This 
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finding does not really come as a surprise, since when expanding the exponential, 
the leading term for det at [D] comes from the term T^'^) in the generalized 
hopping method and from the term Cq in the generalized Fourier method. Thus 
for g > Q the leading terms are missing and both methods cannot be used 
reliably for such values of q. On the other hand, as long as the bound q < Q 
is respected we found that both methods agree well with the results from the 
standard Fourier method. 

5. Properties of canonical fermion determinants 

Having discussed and compared the possible strategies for the evaluation of the 
canonical fermion determinants, we now analyze their properties numerically. 
Except for the discussion of the results for free fermions in the next subsection, 
throughout this section the canonical determinants were evaluated on our pure 
gauge ensembles using the generalized hopping expansion method for quark 
numbers ranging from q = —6 to q = 6. 

5.1. The case of free fermions 

We start our discussion of canonical determinants with the case of free fermions. 
The study of the free case allows for a precise assessment of the three approaches 
discussed in Sections 3.2 - 3.4, and furthermore provides an interesting illus- 
tration of the properties of canonical determinants. In this situation it is rather 
straightforward to evaluate the grand canonical determinant and thus also the 
canonical determinants. For the free case Fourier transformation allows one 
to compute explicitly all eigenvalues of the Dirac operator ([1]) for arbitrary 
boundary conditions. The product of the eigenvalues gives the grand canoni- 
cal determinant and the canonical determinants are obtained easily using the 
Fourier transformation (flOi) . 

It is interesting to analyze how the canonical determinants deti [D]^'^^ are dis- 
tributed as a function of q. For that purpose we consider for one flavor the ratio 
deti[D]('^) / deti[D](''\ i.e., we normalize with respect to the trivial sector q = 0. 
More generally, in the sum for the fugacity expansion the canonical determinants 
are weighted with ex.p{fiqP). In order to directly illustrate the behavior of these 
terms, in Fig. [3] we show the combination exp(/xg/3) deti[L'](''Vdet[£']^°'' for 
/i = and fi = 0.1 (in lattice units). For /x = this reduces to the distribution 
of the canonical determinants. The results we show in Fig. [3] are for lattice size 
8^ X 4 and quark mass m = 0.1 (in lattice units). 

The /i = data in Fig. [3] show that the distribution of the canonical de- 
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Figure 3: Distribution of the (normalized) terms of the fugacity expansion, 
e^''^deti[L>](9)/deti[Li](0) for one flavor of free fermions on an 8^ x 4 lattice 
with mass m = 0.1 (in lattice units). We compare the cases of vanishing 
chemical potential and chemical potential ^ = 0.1 (in lattice units). 

terminants follows a Gaussian type of behavior. The canonical deternninants 
that contribute in the fugacity expansion at /i = are essentially in the sym- 
metrical interval q G [—10,10], since for values of q outside this interval the 
ratio deti[D](''^/ deti[Z)]('') is already very small. As the chemical potential is 
turned on, the canonical determinants receive additional weight factors e'^'^^. 
This leads to an enhancement of the canonical determinants with positive q, 
while the contributions of determinants with negative q are suppressed. The 
whole distribution is shifted and the terms with a sizable contribution to the 
fugacity expansion at fi = 0.1 now come from the interval q G [—7, 17]. 

Changing the value of the mass parameter has only a relatively small effect, 
with the distribution becoming wider as the quark mass in decreased. Increasing 
the spatial volume at fixed (3 leads to a wider distribution. If both the spatial 
volume and the temporal extent (3 are increased at a fixed aspect ratio L//3 and 
fixed mass m, the width of the distribution remains essentially invariant. 



19 



-2e+42 2e+42 -1e+42 



1e+42 -3e+41 



3e+41 



~r 



-4e+38 



T" 



~r 



q = 



T" 



I ' ' ' ' I 

. q = 1 




~r 



q = 2 



1^ 



3e+41 



I ' ' ' ' I 
.Cl = 3 



C>=4 



-- » • 



T — ' — ' — r 
q = 5 



-3e+41 



3e+39 



-3e+39 



-i — , — , — , — |- 
q = 6 



• •• • 



-1e+40 1e+40 -3e+39 3e+39 



4e+38 



Figure 4: Scatter plots in the complex plane for the one-flavor canonical 
fermion determinants deti[Z?]('^) for q = 0,1, 2. ..6 in the center symmetric 
phase. The data are from 1000 configurations on our 8^ x 4 lattices at 
Pg = 7.4. Note that the scale of the plots changes as we increase q. 

5.2. Scatter plots of the canonical determinants 

After illustrating the properties of the canonical determinants in the free case, 
we now switch to our gauge ensembles. We begin our analysis of the canonical 
determinants by showing scatter plots of the values of the one-flavor determi- 
nants deti[L']'^'') for q = 0,1 ...6 using the 8"^ x 4 lattices at two different values 
of the gauge coupling, giving rise to an ensemble in the center symmetric phase 
(/3g = 7.4), and one in the center broken phase = 8.2). 
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Figure 5: Same as Fig. [U but now for an ensemble in the phase with broken 
center symmetry (Pg = 8.2). 

Figure m shows the scatter plots for the ensemble at Pg = 7.4, i.e., in 
the center symmetric phase. For q = the canonical determinant is real (see 
Section 2.2), while for (7 > the values scatter in the complex plane. For 
the two cases with vanishing triality, q = 3 and q = 6, one expects a non- 
vanishing expectation value (deti[D]('^))G according to the discussion of the 
center symmetry in Section 2.4. Thus there must be an anisotropy in the 
distribution in the complex plane. However, such an anisotropy is hardly visible 
for g = 3 and g = 6, which illustrates that computing, e.g., the free energies 
for baryonic states in the confining phase from the corresponding canonical 
determinants is heavily plagued by phase cancellations. 
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In Fig. Owe now show the scatter plots for an ensemble in the phase with 
broken center symmetry {(3g = 8.2). The breaking of the center symmetry 
manifests itself in much more organized scatter plots. Again the q = canonical 
determinant is exactly real, but now the other two cases with vanishing triality, 
q = 3 and q = 6, also have their data points concentrated near the real axis. 
The other cases, q = 1,2,4,5, show the angular distribution familiar from the 
Polyakov loops (compare Fig. [2]), where the data points have phases close to 
the center values 0, 27r/3 and — 27r/3. The pattern is most pronounced for 
small values of q and becomes less clean as q is increased. This behavior 
can be understood due to the fact that the number of windings of the loops 
building up the canonical determinants is given by q: Although the canonical 
determinants always transform according to (i2T]l . loops that wind many times in 
order to give rise to large q amplify the noise introduced by the individual gauge 
configuration, and the pattern reflecting the center transformations becomes 
more blurred with increasing q. Changing the mass has essentially no influence 
on the center pattern. 

It would be very interesting to analyze how the sharpness of the center 
pattern changes as one goes to larger and finer lattices (a question beyond the 
scope of this study). The experience from the Polyakov loop leads one to expect 
that the pattern becomes sharper as the volume is increased and the resolution 
of the lattice becomes finer. Anyway, it is surprising that for the smaller values 
of q the center pattern is so clearly pronounced, much more pronounced than 
for the Polyakov loop as may be seen by comparison with Fig. [2j 

In Section 2.2 we discussed the relation between canonical one-flavor deter- 
minants and the canonical determinants for A'^ > 1 flavors. For completeness in 
Fig. El we show a few scatter plots of canonical determinants in the two-flavor 
case as computed from (fl3|) . In particular we show det2[i?]*^'^'' for g = 0,1,2 
(left to right). The top row is for an ensemble in the center symmetric phase 
{Pg = 7.4), while in the bottom row we use a high temperature ensemble 
(/9g = 8.2). 

It is obvious that for A'^ = 2 qualitatively the behavior is the same as in the 
one-flavor case. In both phases the q = determinants are real. In addition for 
the case of two flavors the values are strictly non-negative as may be seen from 
(fT3]l and (fT4l) . In the center symmetric phase the values for g = 1 and q = 2 
scatter isotropically in the complex plane, while for the center broken phase we 
again find the familiar center pattern. The amplitudes are roughly the squares 
of the amplitudes for the one-flavor case, as is expected from (fTSll where the 
two-flavor case is expressed as a sum of squares of one-flavor determinants. 
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Figure 6: Scatter plot of the canonical two-flavor determinants det2[i?]^''^ 
for q = 0, 1 and q = 2 for our 8^ x 4 lattices. We compare the results for an 
ensemble in the center symmetric phase {(3g = 7.4, top row of plots) to an 
ensemble in the center broken phase {Pg = 8.2, bottom). 



5.3. Implications of the center properties 

In the last section we have illustrated that the center pattern is surprisingly 
strongly pronounced in the high temperature phase. As a consistency check 
whether we understand the connection between the sectors of the Polyakov 
loop as shown in Fig. [2] and the center pattern in Fig. [5l we divide the 1000 
configurations in our Pq = 8.2 ensemble into three subsets according to the 
phase 9p of the Polyakov loop. In each of the three sectors, the real {9p 0), 
and the two complex sectors {9p ±2tt/3), we have roughly a third of the 
total statistics. 

In Fig. [7] we show scatter plots for the canonical q = I (top row) and 
q = 2 (bottom) determinants. In the Ihs. column we display the sector with 
real Polyakov loop, followed by two columns showing the two complex sectors. 

The figure nicely demonstrates, that for g = 1 the center angle follows 
exactly the angles determining the sectors of the Polyakov loop. For q = 2 
the phase of deti[L']*^^^ is complex conjugate compared to the Polyakov loop 
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Figure 7: Scatter plots in the complex plane for the canonical one-flavor 
determinants detifZ)]'-'^^ for q = 1 (top row of plots) and q = 2 (bottom) 
at high temperature. The data are separated with respect to the phase 9p 
of the Polyakov loop. From left to right we show the sectors with 0p ~ 0, 
Op 2n/3 and Op ~ -2n/3. The data are for 8^ x 4 at /3g = 8.2. 

phase because the two-flavor determinant transforms under center rotations with 
= z*. Fig. [7] thus confirms that the center pattern is manifest as expected 
from Eq. (fTTD . 

Let us now discuss an important consequence of the strong center pattern 
observed in the canonical determinants. The transition from the real Polyakov 
loop sector to the sector with Bp ~ 2tt/3 may be implemented by a center 
rotation of the gauge configuration with z = exp(i27r/3) ~ exp(z^p), and 
equivalently for the other complex Polyakov loop sector. According to Eq. (|2T1) 
the canonical determinants transform as deti[D]('^^ z"^ deti[D]^'^^ under this 
rotation. Inserting this into the fugacity expansion (fTT]) we obtain for the grand 
canonical determinant in the different Polyakov loop sectors, 

det[D{n)] = ^zMeti[D](^)e^'?^ . (37) 

Using this equation we now discuss that in the high temperature phase the 
canonical determinant must behave quite differently in the different sectors of 
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6'p ~ 


9p ~ ±27r/3 


0.00 


0.175(15) X lO^y 


0.45(13) X 10^^ 


0.05 


0.212(19) X 10^^ 


0.43(14) X 10^3 


0.10 


0.379(36) X 10^9 


0.37(14) X 10^3 


0.20 


0.447(51) X 10*^° 


0.26(17) X 10^3 


0.40 


0.109(27) X 10^^ 


0.47(47) X 10^^ 


0.60 


0.48(27) X lO'^^ 


0.11(11) X 10^5 



Table 2: Results for (| det[L'(//)]|)G at various values of the chemical poten- 
tial /i in lattice units for our 83 X 4, /3g = 8.2 ensemble. 

the Polyakov loop. Let us first consider the real Polyakov loop sector where 
z = 1. Then the fugacity expansion ([37ll is a sum of essentially real and 
positive terms. The fact that for 9p = the deti[D]^'^'' are essentially real and 
positive follows from the scatter plots in Fig. [5] and [71 

The complex Polyakov loop sectors, on the other hand, are characterized 
by Op ~ 27r/3. As discussed in the last paragraph, these phases can also be 
obtained from the real sector by a center rotation with z = exp(ibi27r/3). 
Inserting these values of z into ([371) leads to relative complex phases in the 
fugacity expansion and thus to cancellations. Consequently one expects that 
the values of the canonical determinants are smaller for the complex Polyakov 
loop sectors. 

In order to study this effect we analyzed the expectation value of the modulus 
of the grand canonical determinant, (| det[D{fi)]\)G, again dividing the gauge 
configurations into the three Polyakov loop sectors. Table [21 shows the corre- 
sponding results for our /3g = 8.2 ensemble. It is obvious from the table, that 
for the complex Polyakov sectors the average size of the grand canonical deter- 
minant is much smaller than for the real Polyakov loop sector. Thus the latter 
sector receives a much larger weight in the path integral and thus is selected 
by the system in the high temperature phase. When analyzing {\det[D{fi)]\)G 
in the low temperature phase, we found essentially no discrepancy between the 
three Polyakov loop sector^. We conclude that in the high temperature phase 
the selection of the real Polyakov loop sector in the dynamical case can be un- 
derstood as a consequence of the center symmetry properties of the canonical 
determinants. 

•^Actually the division of the gauge configurations at low temperatures into three sectors 
according to the phases of the Polyakov loop is somewhat artificial as is obvious from the 
scatter plots for low Pg in Fig. (2] Nevertheless it provides a useful consistency check. 
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5.4. Center symmetry and observables 



Having demonstrated the importance of the center symmetry, we briefly discuss 
its relevance for observables. An observable is a gauge invariant object and as 
such may be expressed as a sum over products of closed loops dressed with 
link variable^, similar to Eq. pTll which gives the loop representation for the 
fermion determinant. Each term in this sum may be classified with respect to 
the total winding number r G Z of its loops. Organizing the loop representation 
of an observable O with respect to the winding numbers, we find 

O = Y^ . (38) 

r 

Since the terms O^*") consist of loops with a total winding number of r, they 
transform under the center transformation (flQll as 



O(r-) , Q{r) _ (gg) 

Let us now analyze the vacuum expectation value of an observable. We find 
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In the second step we inserted the fugacity expansion (fTTll for the determinant 
and the sum p9l) for the observable. The individual terms in the sum transform 
under a center transformation as 

detjv[I?]^'^^0(") ^ zi+'-deiN[Di'^^0^''K (41) 

As long as the center symmetry is unbroken, this implies that the contribution 
vanishes unless 

(g + r)mod3 = , (42) 



which follows from the same reasoning as used in ([22]) . In other words, only 
those terms, where the combined triality of the canonical determinant and the 
observable term O^*") is trivial, give a non-vanishing contribution in vacuum 
expectation values. For high temperature, where the center symmetry becomes 
broken, the three center elements are not summed over in the path integral and 



■^We assume that the fermions are aheady integrated out, such that fermionic observ- 
ables are expressed in terms of quark propagators. 
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the selection rule (j42|) does no longer apply, and O^*") receives contributions 
from all quark sectors in the center broken phase. 

Particularly simple are observables where the sum p9]) reduces to a single 
term with a fixed value of r. Examples are the Polyakov loop (r = 1), the 
Wilson loop (r = 0), the plaquette (r = 0), the dual chiral condensate [llj 
(r = 1) and similar projected observables [12] [r = 1). The whole observable 
in the low temperature phase then receives only contributions from the quark 
sectors q that have {q + r)mod3 = 0. For example the Polyakov loop has 
its leading contribution in the low temperature phase from the sector q = —1, 
while in the high temperature phase the leading contribution comes from q = 
(see also [H [13] for this example). 

5.5. Size distribution properties of canonical determinants 

Analyzing scatter plots in Section 5.3 we established interesting properties for 
the angular distribution of the canonical determinants at high temperature, 
which, as we illustrated, is related to center symmetry. In this section we will 
show that also the modulus of the canonical determinants shows an interesting 
behavior. 

We begin our discussion with Fig. [8l where we plot the average over the 
modulus of the canonical one flavor determinants (| deti[L)]('')/deti[L)]('^)|) as 
a function of q, normalized relative to the q = determinant. In other words 
we study the average size of the canonical determinant for different q. In the 
top panel of Fig. [8] we show for the 8'^ x 4 ensembles the results for the low 
temperature phase, i.e., from /3g = 7.4 up to = 7.8 (which is roughly the 
gauge coupling for the transition). The bottom panel is for high temperature 
{(3g = 7.8, ...8.2). 

The low temperature results presented in the top plot of Fig. [5] show a 
similar behavior as the results for the free case (compare the fi = data in 
Fig. E]). The distribution is symmetrical around g = (as it must be since 
deti[Z)](~'') = (deti[D] '■''))*) and has a Gaussian-like distribution as a function 
of q. The distribution is rather narrow for the lowest temperature (/3g = 7.4) 
and widens as (3g and thus the temperature are increased. 

As a consistency check we summed up the fugacity expansion (fTTjl with 
the available terms, q G [—6,6], and compared it to the exact result for the 
grand canonical determinant. For the lowest values of Pg, where, as Fig. [8] 
shows, the canonical determinants become small very quickly, we found excellent 
agreement, and the fugacity expansion with q G [—6,6] reproduces the grand 
canonical ^ = determinant with an error smaller than 10^^ (on 8'^ x 4). When 
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Figure 8: The average size of the determinant normahzed relative to the 
q = case, (| deti[Z)]('^)/ deti[L']('') |), as a function of q. In the top plot we 
show the results for the gauge couplings (3g = 7.4, 7.5, 7.6, 7.7, 7.8 (T < Tc), 
while the bottom is for (3g = 7.8,7.9,8.0,8.1,8.2 (T > Tc). The data are 
from our 8^ x 4 ensembles. 



increasing fi, which shifts the distribution towards larger q (compare Fig. [3]), we 
found that the error remains below 1% up to ^ = 0.35 (in lattice units). This 
finding is an interesting result which shows the fast convergence of the fugacity 
expansion at low temperatures and moderate lattice volumes. Beyond that 
it is also an important consistency and accuracy check for our evaluation of 
the canonical determinants. We stress that the achieved accuracy was not only 
established for the ensemble average, but holds also for individual configurations. 

Let us now come to the bottom plot in Fig. [51 which shows our results 
in the high temperature phase. It is obvious that the trend towards a wider 
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distribution at increasing temperature (increasing fJc) continues, although for 
the largest values of fSc the widening slows down. Most interesting, however, 
is the observation that the shape of the distribution is no longer Gaussian, but 
a relative enhancement of the determinant values for even q becomes visible for 
the highest values of Pa- 
in order to study this enhancement of even quark numbers further, in Fig. [9] 
we compare the distribution for the low and high temperature phases now using 
the two ensembles on the larger lO'^ x 4 lattices. While the canonical deter- 
minants in the low temperature phase [Pg = 7.5, T/Tc ~ 0.77) again show a 
Gaussian-like distribution, for the high temperature phase {Pg = 8.1,T/Tc ~ 
1.32) the enhancement of the even quark numbers is much more pronounced for 
the larger 10^ x 4 volumes shown in this figure. On the other hand, we found 
that on our smaller 6^^ x 4 volumes the effect vanishes. Thus the enhancement 
of even quark numbers is an effect that needs a minimal spatial volume, and 
probably would become more pronounced if the volume was increased further. 

The relative enhancement of the canonical determinants for even quark 
numbers may be interpreted as an enhancement of meson states, of diquark 
states, et cetera in the high temperature phase, while single quark states, baryon 
states, et cetera are suppressecG. A possible formation of diquark correlations 

^We stress at this point, that of course there are particle anti-particle creation and 
annihilation processes in each sector. The terms "quark sectors", "diquark states" et 
cetera, thus always refer to the net surplus of these objects. 
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at not too large temperatures has been suggested [14] and preliminary lattice 
studies may be found in the literature [TSj. Although it certainly needs to 
be improved (larger volumes, more systematical scanning of temperatures and 
masses), our study provides further ab-initio support for such a picture and 
generalizes it to all even quark number sectors. Of course a crucial question, 
which will have to be addressed in future work, is whether the effect persists 
when a better sampling beyond the pure gauge ensembles is being used. Such 
an analysis, as well as the question what happens at larger values of T/Tc (here 
we have T/Tc ~ 1.32) is planned for the future. 

We stress again at this point, that in the low temperature phase also a 
suppression of certain quark sectors is manifest - the sectors with non-trivial 
triality. At low temperature, however, the mechanism is a different one, as we 
discussed in Section 2.4: The unbroken center symmetry leads to the survival 
of only the sectors with vanishing triality, i.e., mesons, baryons et cetera. This 
happens through the averaging over the phases of the canonical determinants. 
Since in our observable (| deti deti[-D]('^)|) the phases are neglected, the 
suppression of the sectors with non-vanishing triality is not manifest in the 
corresponding low temperature data in Figures[8]and [9l 

6. Summary and outlook 

In this article we discuss the numerical evaluation and properties of canonical 
fermion determinants using Wilson fermions on pure SU(3) configurations. We 
obtain a considerable speedup in the evaluation of the canonical determinants by 
implementing a domain decomposition approach which leads to a dimensional 
reduction of the fermion determinant [6J. 

Based on this dimensional reduction we compare three different strategies for 
evaluating the canonical determinants: The conventional Fourier transformation 
approach, a generalized hopping expansion, and a generalized Fourier method, 
which combines a Fourier transform of the logarithm of the determinant with 
expansion techniques. Our comparison shows that all three methods work, 
are sufficiently accurate for the volumes, quark numbers q and couplings we 
considered, and allow for a systematical improvement of the accuracy if needed. 

We evaluate the canonical determinants with quark numbers q G [—6,6], 
which for our lower temperatures and not too large chemical potential is suf- 
ficient for a very accurate representation of the grand canonical determinant 
through the corresponding fugacity expansion. The obtained accuracy of the 
fugacity expansion is an important test for the correctness of our implementa- 
tion. For larger chemical potential and high temperatures more values of the 
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quark number q would be needed. 

Subsequently we explore the properties of the canonical determinants. Scat- 
ter plots of the determinant values show that in the low temperature phase the 
determinant values are distributed essentially isotropically in the complex plane. 
This holds also for the sectors with vanishing triality (except for the zero quark 
sector), illustrating that the evaluation of the free energy of, e.g., a baryon from 
the corresponding q = 3 canonical determinant is heavily plagued by phase 
cancellations. 

In the high temperature phase we observe a very pronounced center pattern 
of the canonical determinants and we show that this pattern behaves exactly 
as expected from the center transformation properties of the canonical deter- 
minants. 

Based on the observation that the center symmetry is so strongly mani- 
fest in the canonical determinants, we explain and numerically confirm that the 
grand canonical determinant is much smaller for high temperature gauge con- 
figurations with complex phases of the Polyakov loop. Thus the fact that the 
Polyakov loop comes out real in dynamical simulations may be understood as 
a consequence of center symmetry and its breaking. The discussion of the role 
of the center symmetry is furthermore extended to observables, giving rise to 
selection rules for the relevant quark sectors in the low temperature phase. 

Analyzing the distribution of the canonical determinants as a function of the 
quark number, we find for low temperatures a Gaussian-like distribution which 
widens as the temperature is increased. At high temperatures an interesting ef- 
fect is observed: Canonical fermion determinants with even quark numbers q are 
enhanced in the high temperature phase. This might indicate an enhancement 
of diquark states et cetera. 

We conclude with remarking again that an important caveat must be kept 
in mind: The numerical results were obtained using ensembles of pure gauge 
theory. Although this is formally correct, in a numerical simulation it could be 
the case that the distribution of the generated gauge configurations has only a 
small overlap with the main contributions to the canonical determinants, and 
very large ensembles might be needed for accurate results. This is an open issue 
which might be answered by using dynamical ensembles and suitable reweighting 
- a question that we plan to address in the future. 
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